Quantum properties of transverse pattern formation in second-harmonic generation 
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We investigate the spatial quantum noise properties of the one dimensional transverse pattern 
formation instability in intra-cavity second-harmonic generation. The Q representation of a quasi- 
probability distribution is implemented in terms of nonlinear stochastic Langevin equations. We 
study these equations through extensive numerical simulations and analytically in the linearized 
limit. Our study, made below and above the threshold of pattern formation, is guided by a micro- 
scopic scheme of photon interaction underlying pattern formation in second-harmonic generation. 
Close to the threshold for pattern formation, beams with opposite direction of the off-axis criti- 
cal wave numbers are shown to be highly correlated. This is observed for the fundamental field, 
for the second harmonic field and also for the cross-correlation between the two fields. Nonlinear 
correlations involving the homogeneous transverse wave number, which are not identified in a lin- 
earized analysis, are also described. The intensity differences between opposite points of the far 
fields are shown to exhibit sub-Poissonian statistics, revealing the quantum nature of the correla- 
tions. We observe twin beam correlations in both the fundamental and second-harmonic fields, and 
also nonclassical correlations between them. 



I. INTRODUCTION 

Pattern formation has been an active area of research 
in many diverse systems Q. Numerous similarities to 
pattern formation in other systems have been reported in 
recent studies in nonlinear optics §~q|- However similar, 
nonlinear optics also displays properties that are wholly 
unique due to the relevance of quantum aspects in optical 
systems, one manifestation of this is the inevitable quan- 
tum fluctuations of light. In the last decade an effort has 
been made to study the interplay in the spatial domain 
between optical pattern formation, known from classical 
nonlinear optics, and the quantum fluctuations of light 
New nonclassical effects such as quantum entan- 
glement and squeezing in patterns were predicted |&^|. 
Another interesting example is the phenomenon of quan- 
tum images: below the instability threshold, information 
about the pattern is encoded in the way the quantum 
fluctuations of the fields are spatially correlated Il0| . 

Nonlinear x' 2 '-materials immersed in a cavity have 
shown most promising quantum effects. A paradigm 
of spatiotemporal quantum behavior has been the opti- 
cal parametric oscillator (OPO), which despite its strik- 
ing simplicity is able to display highly complex behav- 
ior In the degenerate OPO, pump photons 
are down-converted to signal photons at half the fre- 
quency and with a high degree of quantum correlation. 
This might be attributed to the fact that the signal 
photons are created simultaneously conserving energy 
and momentum, leading to the notion of twin photons. 
In the opposite process of second-harmonic generation 
(SHG) fundamental photons are up-converted to second- 
harmonic photons at the double frequency. On a clas- 



sical level, both the OPO and intra-cavity SHG display 
similar spatiotemporal behavior. The essential difference 
between them is that in the OPO an oscillation threshold 
for the process exists, which simultaneously acts as the 
threshold for pattern formation. On the contrary, SHG 
always takes place no matter the strength of the pump 
field, but there is a threshold that marks the onset of pat- 
tern formation. This gives pronounced differences with 
the OPO in the linearized behavior below the threshold 
for pattern formation. In the OPO the pump and the 
signal fields effectively decouple and only the latter be- 
comes unstable at threshold. At a microscopic level, the 
behavior of the OPO close to the threshold can be un- 
derstood in terms of a unique process in which a pump 
photon decays into two signal photons with opposite wave 
numbers. In SHG the fundamental and second-harmonic 
fields are coupled and both become unstable at threshold. 
This complicates the picture mainly by the number of mi- 
croscopic mechanisms that are relevant to describe the 
pattern formation process. But this complexity, on the 
other hand, is likely to generate interesting correlations 
between the fundamental and the second-harmonic field. 
Recently, transverse quantum properties in the singly res- 
onant SHG setup were investigated |Q. There, squeez- 
ing in the fundamental output was observed close to the 
critical wave number, but since the second-harmonic is 
not resonated the question of possible correlations be- 
tween the two fields was not addressed. However, since 
the second harmonic in the singly resonant case is given 
directly as a function of the fundamental, correlations 
similar to the ones observed in the fundamental should 
be expected. In this paper, we will consider the case of 
doubly resonant SHG with the aim of investigating the 
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spatial correlations not only within each field (fundamen- 
tal field and second- harmonic field), but also between the 
two fields. 

For this purpose we use the formalism of quasi- 
probability distributions jl5fl . Choosing the use of the 
Q representation we are able to derive a set of nonlin- 
ear Langevin equations that describes the time evolution 
of the quantum fields in the SHG setup (Sec. |J). -"- n 
Sec. ID, the linear stability analysis of this system will 
be discussed and a proper regime of parameters specified, 
for which the formalism adopted here is applicable. Sec- 
tion |y] will be devoted to an analysis, on a microscopic 
level, of the implications of the three-wave interactions in 
the nonlinear crystal. These considerations allow to iden- 
tify the most important spatial correlations expected in 
this two-field system and to define suitable quantities to 
be calculated. In particular, we will focus on equal time 
correlation functions of intensity fluctuations and we will 
study photon number variances when looking for non- 
classical features of the intra-cavity fields. A systematic 
study of the spatial correlations is presented first through 
analytical results in the framework of a linearized theory 
below the threshold for pattern formation (Sec. [v|), and 
also through extensive numerical simulations of the non- 
linear Lang evin equations reported below (Sec. [Vl|) and 
above (Sec. VII) the t hreshold for pattern formation. We 
conclude in Sec. VIII. 
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FIG. 1. The model setup in a top-view. 



II. NONLINEAR QUANTUM MODEL FOR 
INTRA-CAVITY SHG 

We consider a nonlinear \^ -material with type I 
phase matching immersed in a cavity with a high reflec- 
tion input mirror M\ and a fully reflecting mirror M 2 at 
the other end, cf. Fig. [I]. The cavity is pumped at the 
frequency w\ and through the nonlinear interaction in 
the crystal photons of frequency uj 2 — 2a>i are generated. 
This is the process of SHG. The cavity supports a discrete 
number of longitudinal modes, and we will consider the 
case where only two of these modes are relevant, namely 
the mode tt>i, cav closest to the fundamental frequency and 
W2,cav closest to the second-harmonic frequency. In the 
setup shown in Fig. [j] ui 2 . ca ,v = 2wi, cav , but we will allow 
the cavity resonances to be independent in order to con- 
trol the detunings individually. The pump beam propa- 
gates along the z-direction and using the mean field ap- 
proximation, variations in the z-direction are averaged 
out. This approach is justified as long as the losses and 
detunings are small. Due to diffraction the transverse 



section perpendicular to the z-direction spanned by the 
xy-plane also comes into play. We consider the simple ID 
case where only one of the transverse directions is rele- 
vant, so variations along the y-direction are neglected 
and only the x-direction is taken into account (this could 
easily be achieved experimentally by using a crystal with 
a small height). Let A\{x,t) and A 2 (x,t) denote the 
ID intra-cavity boson operators of the fundamental field 
(FH) and second- harminic field (SH) , respectively. They 
obey the following equal time commutation relation 

[A i (x,t),A](x',t)]=S ij S(x-x'), i,j = 1,2. (1) 

The Hamiltonian operator describing SHG including 
diffraction can be written as done in Ref. Jl^] for the 
OPO, 

H = if tree + H int + H CKt , (2) 
where the free Hamiltonian is given by 

f ( c 2 d 2 \ - 

fffree = %] dxA[(x,t) f - __J A 1 { X ,t) 

+hJdxAl(x,t)(-S 2 -—-^jA 2 (x,t). (3) 

Here Sj = ujj — Wj. cav are the detunings from the near- 
est cavity resonances, d 2 /dx 2 describes the diffraction, 
and c is the speed of light. The interaction Hamiltonian 
describes the nonlinear interaction in the material 

Hint = ^ J dx (A 2 (x,t)(Al(x,t)) 2 - H.c.) , (4) 

where g is the nonlinear coupling parameter proportional 
to the x^- n onlinearity f the crystal. The external 
Hamiltonian describes the effects of the pump injected 
into the cavity at the fundamental frequency, which is 
taken to be a classical quantity £;„, so we have 



.Hext = ih / dx (& in A\{x,t) - £* a Ax(x,t) 



(5) 



Then the master equation for the density matrix p in the 
interaction picture is given by 



f = -i[^] + (X 1+ L 2) p. 



(6) 



The cavity losses are assumed to occur only through 
the input coupling mirror to the external continuum of 
modes, and are here included through the Liouvillian 
terms 



L 



,p = j dx^ j (2A j (x,t)pAt(x,t) 

-pAX(x,t)Aj - At( x ,t)A j (x,t)p), 



(7) 



where 7j are the cavity loss rates. Here we have assumed 
that thermal fluctuations in the system can be neglected. 
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Using the standard approach of expanding the den- 
sity matrix into coherent states weighted by a quasi- 
probability distribution function, the master equation (||) 
is mapped onto a functional equation, depending on the 
order for creation and destruction operators M,fL7j. For 
a Hamiltonian that is quadratic in the field operators 
this results in a Fokker-Planck equation, implying that 
the dynamical evolution of the distribution function may 
be modelled by an equivalent set of classical stochastic 
Langevin equations. However, due to the contributions 
of higher order to the Hamiltonian (|4|) problems may 
arise. When using the Wigner representation the evo- 
lution equation of the quasi-probability functional con- 
tains third order derivatives, which means that no equiv- 
alent Langevin equations can be found. These third order 
terms have been shown to model quantum jump processes 
|l8l| . In some cases these terms can be neglected lead- 
ing to an approximate Fokker-Planck form. When using 
the P or Q representations problems of negative diffu- 
sion in the Fokker-Planck equation come into play Jl5| . 
To avoid negative diffusion in the P representation, some 
techniques have been developed where the phase space 
is doubled jl9|,|o|], but then numerical problems due to 
divergent stochastic trajectories generally appear pl"| , f22[ . 
We choose here to use the Q representation which in a 
restricted domain of parameters has a nonnegative diffu- 
sion matrix and has been shown to be a useful alternative 
in the similar problem of calculating nonlinear quantum 
correlations in the OPO p3|] . The Q representation has 
no singularity problems, is bounded and always nonneg- 
ative. 

Introducing ctj and a* as the c-number equivalents of 
the intra-cavity boson operators Ai and Aj , the evolution 
equation for the quasi-probability distribution function 
Q(a) is 



gate* -i^JL-l^ 




■x \ _l # 2 . c 2 d 2 g d 2 
^[(72 - + -a, - ^] - -a 2 ^ 

+ ^^ + ^0a^ 2 +C - C ) Q ^ (8) 

with a = {ax, a*, a 2 , cn 2 }. This is just an extension to 
the diffractive case of the result obtained by Savage |24j . 
Equation (||) has the form of a Fokker-Planck equation, 
and it has positive diffusion if 



M < 2 



7l 



(9) 



As shown below, it is possible to fix the parameters of the 
system in such a way that the stable solution for the SH 
field is well below the value 2^/i/g. Fluctuations around 
this stable solution are small, so that the probability vi- 
olating the condition (^|) is almost zero. Neglecting then 
stochastic trajectories violating this condition, we may 



write a set of equivalent Langevin stochastic equations by 
applying the Ito formalism for the stochastic integration 
p5| . We then obtain the following nonlinear Langevin 
equations 

d t ai(x, t) = (—71 + iS 1 )ai(x, t) + ga\{x, t)a 2 (x, t) 
c 2 d 2 



+*o — ^~2 a xOM) +£m + V 2 7i£i( : M), (10a) 
d t a 2 {x, t) = (-72 + iS 2 )a 2 (x, t) - ^a 2 (x, t) 

+ C a 2 <X2{x,t) + y/*nb{x,t), (10b) 

with multiplicative Gaussian white noise sources corre- 
lated as follows 



(e*(z, t)^( x ', t')) = s^x - x ')s(t - ^ 

ga 2 (x,t) 



(11a) 
(lib) 



(^(x,t)^(x',t')) 
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8{x - x')S(t - t'). (11c) 



We rescale space and time according to 

i=tji, x = x/l d , (12) 
where Id is the characteristic length scale given by 



7 2 - 
l d — 



271^1 

We also normalize the fields and noise according to 



Aj(x,t) = aj(x,t) — , ij(x,t) = £j(x,t) 4/ — , 

71 V 71 



(13) 



E = £ u 



7?' 



(14) 



This allows to rewrite the Langevin equations in dimen- 
sionless form: 

diA x {x,l) = (-1 + i&i)Ai{x,T) + Al(x,t)A 2 (x,t) 
dx* V n t h 



+ i l^ 2 Mx,t) + E + x j— &(£,*), (15a) 



diA 2 {x,t) = {-j + iA 2 )A 2 {x,t) - \a\{x,1) 

id 2 , ,_ ~. / 27 - . ~. 
+ 2gtf -,Mx,t) + J^ 2 (x,t), 



(15b) 



where 7 = 72/71 and Aj — Sj/ji, and E may be taken 
real. Moreover we have introduced 



rhh 



l 2 l d 



(16) 



which in the OPO coincides with the number of photons 
in the characteristic "area" Id required to trigger the os- 

— 1/2 

dilation. The noise strength is seen to scale like n th 
The normalized noise sources are correlated by 
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(£(£,*)&(*'>*')) = W - x')S(t-F) 
(&(£,*)&(£',?)> = 0, 

A 2 (5,t) 



«i(5 > t)fi(5',f)> = 



(17a) 
(17b) 

<y(x-5')*(*-?)- ( 17c ) 



These noise sources turn out to only to be defined for 

\A 2 (x,t)\ < 2, (18) 

which coincides with the condition (^) for a positive dif- 
fusion expressed in terms of the rescaled fields. 

In the following the tildes are dropped, and only nor- 
malized dimensionless equations are considered. We will 
also use the terminology u> = a>i and 2u> = lu 2 . 



III. LINEARIZED EQUATIONS AND 
BIFURCATION DIAGRAM 

In this section we consider the linearization of the 
nonlinear Langevin equations in the Q representation 
around the homogeneous steady state solutions below 
the threshold for pattern formation. This approach re- 
lies on the assumption that the fluctuations are small 
with respect to the field mean values, and therefore we 
expect this approach to break down close to the insta- 
bility threshold. We will come back later (Sec. VI B) to 
the question of the validity of the linear approximation. 



We write the fields as Aj(x,t) = Aj + [3j(x,t), where 
flj(x,t) represent the fluctuations around Aj. The clas- 
sical homogeneous values Aj of the fields given by the 
homogeneous steady state solutions of the deterministic 
limit (nth 
this in Eqs. 
equations 



of Eqs. (|15|), as found in Ref. ||. Using 
we find the following set of linearized 



dtfo{x,t) = (~l + iA 1 )/3 1 {x,t)+A 2 f3* 1 (x,t) 



+ Atfc(x,t) 



d 2 



We would like to mention that the Wigner representa- 
tion, in the linear regime, would lead to equivalent results 
without suffering from any limitation since it satisfies a 
Fokker-Planck equation for any value of \A 2 \. However, 
for the sake of a consistent presentation of our results we 
have chosen to consider the Q representaion also in the 
linear case. 

It is instructive to introduce the spatial Fourier trans- 
form of the fluctuations 



d.v 



0j(x,t)e ikx , 



(22) 



which physically represents the amplitude of the fluctu- 
ations in the far field. Considering Eqs. ([l9]) and their 
complex conjugates, it is readily shown that these am- 
plitudes (3j(k,t) fulfill a set of equations which can be 
written in the following matrix form 



M(fc) 



M(fc) 



/ ft (M) 

0a(M) 

V/3 2 *(-M) 

y/r&(-k,t) 

<Tl(fc) 

A* 2 
-At 




/ ft (M) 

/3a(M) 

\/?f(-M) 




(23a) 



(23b) 



where a\(k) = — l + i(Ai — k 2 ) and a 2 (k) = —-f + i(A 2 — 
k 2 /2) have been introduced and each noise term rjj(k,t) 
is the Fourier transform of the noise term appearing in 
the real space linearized Langevin equations (|l9|). Their 
correlations are given by 



£iOM), 

n t h 



(19a) 



(n*(k, t) Vj (k', 0) = M(* - k'Wt - 
A 2 



d t (3 2 (x,t) = (-j + iA 2 )0 2 (x,t)-Aif3 1 (x,t) 



i d 2 



nth 



(19b) 



( m (k,t) m (k',t')) 
(n 2 (k,t) m (k',t')) = 0. 



5(k + k')S(t-t'), 



(24a) 
(24b) 
(24c) 



The correlations of the stochastic sources £i(x, t) in the 
linearized limit become 

(Q(x,t)^(x',t r )) =6 lJ 5(x-x')5(t~t r ), (20a) 

(^(x,m 1 (x',t'))=-^S(x-x')5(t-t'), (20b) 

(6(x,t)6(x',t')) =0. (20c) 

With A 2 being merely a constant, the noise in the linear 
approximation is not multiplicative any more. However, 
as in the nonlinear equations we have the restriction 



\A 2 \<2. 



(21) 



The linear stability of the classical equations obtained 
as the nth — > oo limit of Eqs. ( |l9| ) was investigated by Et- 
rich et al. H . A rich variety of instabilities was shown to 
exist: A self-pulsing instability, that leads to oscillations 
of the homogeneous steady states without any transverse 
structure, was present for all parameters. The oscilla- 
tory transverse instability leading to patterns traveling 
in space and time was only present for certain parame- 
ters and branched out from the self-pulsing instability. 
Bistability was demonstrated for large detunings of same 
sign and for 7 small. Most importantly, for all param- 
eters also stationary transverse instabilities were found 
to exist, i.e. instabilities at a critical transverse wave 
number k = k c and with zero imaginary eigenvalue. It 
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was shown that stripe-type solutions exist but are always 
unstable, and numerical simulations showed that instead 
hexagons are the dominating stationary transverse insta- 
bility. The ID configuration we have chosen to consider 
here has the advantage that the pattern always will be 
a stripe and therefore leads to simpler interpretation of 
the correlations. We will choose a range of parameters in 
which the stationary transverse instability is accessible as 
the primary bifurcation. This bifurcation is supercritical 
in the ID model. 
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FIG. 2. Stability diagram for Ai = 2.0 and 7 = 0.5, show- 
ing transverse stationary instability (solid line), transverse os- 
cillatory instability (dashed line) and self-pulsing instability 
(dotted line). 
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FIG. 3. Transverse instability for Ai = 2.0 and 7 = 0.5 
shown for the intra-cavity second-harmonic field, along with 
the limit for the Q representation, \Ai\ < 2. 

The choice of parameters must take into account the 
requirement of applicability of the Q representation. One 
finds that Eq. ( pT| ) can only be satisfied for Ai > pfjj . 
Using the expressions presented in Ref. |(| and fixing 
Ai = 2.0 and 7 = 0.5 we obtain the bifurcation dia- 
gram shown in Fig. || [^tJ. We observe that for A2 < 
it is possible to obtain stationary patterns (solid line) as 
the primary bifurcation at a critical value of the pump, 



E t ; increasing the pump beyond E t eventually the sys- 
tem will also become self-pulsing unstable (dotted line). 
For A2 > the transverse oscillatory bifurcation (dashed 
line) is the primary one, and therefore travelling waves 
are observed in this region. The bistable area is located 
for A2 > 8.3 and hence beyond the range shown here. 

Expressing the onset of transverse instability, seen in 
Fig. |2j, in terms of the intra-cavity value of the SH we 
have the bifurcation diagram for the transverse instabil- 
ity shown in Fig. ||. We see that for A2 < we are well 
below the limit for positive diffusion (pl[). Therefore, the 
probability of trajectories violating the condition (|l8| ) of 
the nonlinear equations is almost zero. For A2 > 0, in- 
creasing 7 or decreasing Ai towards zero, this threshold 
gets closer to I.A2I = 2. 

We will therefore use the parameters Ai = 2.0, A2 = 
—2.0 and 7 = 0.5 in the rest of this paper, which gives 
a pattern formation threshold of E t = 7.481757 and a 
critical wave number k c = 1.833. The noise strength is 
set to n t h = 10 8 which is a typical value for the cavity 
setup discussed here |2q] . 
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FIG. 4. Numerical simulation of the Langevin equations 

above threshold with E/E t = 1.01 and L = 102.84. Left: 

The absolute value of the near field of the FH (above) and SH 

(below). Right: Far field average intensity of FH, (\Ai{k)\ 2 ). 

The far field of the SH shows a similar structure. 

The main task of the following section is to identify the 
most important correlations we expect to find in the sys- 
tem. For this purpose it is useful to have a good knowl- 
edge of the spatial structures that emerge in the system. 

Numerical simulations |2j| of the nonlinear Eqs. ( |T5| ) 
confirmed the instability at a finite transverse wave num- 
ber k — k c predicted by the linear stability analysis. 
Above the threshold for pattern formation modulations 
was observed around the steady state with wavelengths 
corresponding to k c . This is shown in Fig. |i] where the far 
field intensity shows distinct peaks at k = 0, correspond- 
ing to the homogeneous background, and at k = ±fc c 
corresponding to the modulations observed in the near 
field, as well as higher harmonics. 

Below threshold the quantum noise will excite the least 
damped modes and precursors of the spatial pattern are 
observed. This is shown in Fig. ^| where a space-time 
plot is presented for the FH near and far field. Clearly a 
stripe-type pattern is formed, but as time progresses the 
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noise diffuses the pattern so that averaging over 

time will wash out this emerging structure and a spatially 
homogeneous near field will remain. On the contrary, as 
we will show the spatial correlation functions do encode 
precise information about the emerging pattern, even af- 
ter this time averaging has been carried out, as illustrated 
through the concept of quantum images pQ| . 
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FIG. 5. Numerical simulation with E/E t = 0.9999 and 
L = 103.057, showing the space-time evolution of \A\\ in the 
near field (left) and far field (right). A similar behavior is 
seen for the SH. 



IV. CORRELATIONS, PHOTON INTERACTION 
AND PATTERN FORMATION 

Our general objective is the investigation of the spatial 
intra-cavity field correlations emerging in this system as 
a result of the coupling of FH and SH field through the 
nonlinearity of the crystal, and the implications of the 
spatial instability on these correlations. This study has 
a two-fold purpose: First, to obtain a precise picture on 
how pattern formation occurs in cavity SHG. In particu- 
lar, we will aim at identifying the relevant mechanisms, 
in terms of elementary three- wave processes that are im- 
portant for the understanding of the intra-cavity field 
dynamics. Secondly, it will be interesting to investigate 
whether these correlations are the manifestation of non- 
classical states of the fields. Such states are identified by 
investigating the statistics of the intra-cavity intensities, 
looking in particular for possible sub-Poissonian features 

A. Photon interaction 

We will start by investigating the equal time correla- 
tions between intensity fluctuations at different points in 
the far field. The intensity of each field being directly pro- 
portional to the number of photons in the corresponding 
mode, we can relate the intensity fluctuations to the cre- 
ation or destruction of photons. The idea is that the way 
these fluctuations are correlated gives information about 
the microscopic mechanisms that take place in the cavity 
and, ultimately, that are involved in the pattern forma- 
tion process. Generally speaking, a positive correlation 
tells us that there should exist a coherent mechanism that 



creates simultaneously the corresponding photons. The 
following normalized correlations are considered 

mtMhm , (25) 

where the superscript n denotes normalization. The in- 
tensity fluctuations are given by 6Nj(k,t) = Nj(k,t) — 
(Nj(k,t)), which involves the photon number operator 
Nj(k,t) = Ajj(k,t)Aj(k, t). The brackets denote quan- 
tum mechanical averages (expectations) of the operators, 
which in our approach based on numerical simulations of 
equivalent c-numbers will be translated into an average 
over time. The normalization of the correlations implies 
that C"j(k,k r ) — 1 for perfectly correlated fluctuations, 
whereas C™- (&,£;') = — 1 will be the signature of perfect 
anti-correlation between the intensity fluctuations. As 
usual, the absence of any correlation will translate into 
a vanishing correlation function C^(k,k') — 0. In the 
following we will refer to Cj\(fc, k') and C^ik, k') as self- 
correlations (between different modes of a given field) 
and to C" 2 (fc, k') as cross-correlations (between modes in 
different fields). 

As a guideline for the investigation of the properties 
of these correlation functions, the first step consists in 
identifying the basic photon processes when the system 
is taken close to a transverse instability. These pho- 
ton processes must obey the standard energy and mo- 
mentum conservation laws. Whereas the former merely 
implies that each elementary process must connect one 
SH photon with two FH photons, the latter will trans- 
late into a condition on the transverse wave numbers. 
Keeping in mind that the cavity is pumped with an 
homogeneous field at the frequency to, the first pro- 
cess to consider consists in two homogeneous FH pho- 
tons, [u>](k = 0) = M(0) combining to give one ho- 
mogeneous SH photon, [2w](0), what will be written as 
M(0) + M(0) -> [2w](0). This is encoded in the Hamil- 
tonian term Af A\ in Eq. (§. The inverse process, which 
corresponds to the degenerate OPO process, also takes 
place in the system, as shown by the presence of the term 
(A{) 2 j42 in Eq. (Q). Elaborating on these considerations 
we propose the scheme in Fig. ^| as the simplest way of 
obtaining a pattern in both fields. 

1) The first step is the basic SHG channel where two 
homogeneous FH photons give a SH photon and 
vice versa, i.e. the channel [w](0) + M(0) <-> 
[2w](0). It is important to realize that fluctuations 
around the steady state are considered, hence it is 
not considered how the FH photons combine to give 
the steady state SH photons via the channel above, 
but rather how the fluctuations invoke the channel 
beyond this. 

2) The second step is the down-conversion of a SH 
photon into two FH photons. Momentum conser- 
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vation in the process implies that the two FH pho- 
tons have the same value of the transverse wave 
number but with opposite signs. These are called 
twin photons since an emission of a [ct>] (+£;') photon 
must be accompanied by an emission of a k') 
photon, and they therefore show a high degree of 
correlation. This channel written as [2cj](0) <-> 
[u](—k') + [oj](+k') generates off-axis FH photons. 

3) Off-axis SH photons are obtained by combining 
the created off-axis FH photon from step 2) with 
a photon from the homogeneous background to 
give a SH photon, which by momentum conserva- 
tion must have the same wave number as the off- 
axis FH photon. This channel can be written as 
[w](0) + [w] (+*:') «-» [2w] (+*')• 



3) 



1) 



o(0) s 

- - -> _2oj(0)_ x 



©(0)' 




2©(+k') 



2©(-k') 



FIG. 6. The basic picture of pattern formation on a micro- 
scopic level through SHG. The single arrows (— >) symbolize 
FH photons, while double arrows (=>) symbolize SH photons. 
The dashed arrows are photons from the homogeneous back- 
ground. 



Of course, these are not the only three- wave processes 
which are kinematically allowed in the nonlinear crystal, 
since the interaction Hamiltonian (ji|) induces any pro- 
cess of the form M(fc') + [tu}(k") [2uj}(k' + k"), with 
arbitrary wave numbers kl and k". In fact, the basic 
scheme we propose in Fig. ||only takes into account those 
three-wave processes which involve at least one photon 
of the homogeneous background fields. Empirically, this 
choice is motivated by the observation that below the 
threshold these are the only field modes that are macro- 
scopically populated, so that any process involving them 
should be stimulated in analogy to what occurs in stan- 
dard stimulated emission. Formally, the selection of these 
particular elementary processes corresponds precisely to 



the approximation made by linearizing the field equa- 
tions around the steady state solution. As can easily be 
checked, the full equations for the far field fluctuations 
contain additional terms quadratic in the fluctuation am- 
plitudes, which indeed account for othe r thr ee- wave pro- 
cesses. Linearizing we are left with Eq. (23a), which only 
take into account the processes represented by step 2) 
and step 3). These processes translate into nondiagonal 
elements of the matrix M(fc) of the linear system, and 
as a consequence, for any value of k, the time evolution 
of the four amplitudes fli(k,t), 0i(—k,t), 02(k,t) and 
02(—k,t) will be coupled. This coupling is expected to 
translate into correlations between the intensity fluctua- 
tions 6h(k), 6h(-k), 5h(k), 51 2 {-k). 

This preliminary observation already allows to give a 
more explicit interpretation of the basic scheme of Fig. |(]. 
Splitting the dynamics of the intra-cavity fields into in- 
dependent elementary steps, as suggested in the discus- 
sion of Fig. H, would not explain any correlations either 
between [w](fc') and [2cj](— k') nor between [2u>](fc') and 
[2a;] (—&')• Hence the inspection of the linearized equa- 
tions shows that the interpretation of Fig. |^ in terms of a 
cascade is too naive. Instead, we have to understand step 
2) and 3) as two coherent, joint processes, which generate 
simultaneously correlations between the 4 modes [o>] (£;'), 
[u)](—k r ), [2u](k') and [2u](—k'). Finally, it is important 
to stress that the linearized analysis does not predict any 
correlation between intensity fluctuations in field modes 
with wave numbers of different modulus. Mathemati- 
cally, this is due to the fact that in the linear approxima- 
tion all correlation functions (Eq) have the structure 



Cij(k, k') = C\;\k)5 2 {k - k') + C { J\k)5 2 {k + k'), 



(26) 



as will be shown in the next section. Close enough to 
threshold, however, this will not be true any more be- 
cause of the emergence of additional correlations of non- 
linear nature. 

Let us finally briefly address the fundamental differ- 
ence between OPO and SHG: Whereas in SHG, the two 
fields A\ and A2 are always nonzero regardless of the 
pump level, in the OPO case below the oscillation thresh- 
old A2 is fixed by the pump and Ai — 0. Considering 
the scheme presented in Fig. [| the vanishing of A\ im- 
plies that there is no macroscopic population of the mode 
[w](0) and therefore, step 3) of Fig. || is not present. 
The route to pattern formation simply consists of step 
2) in Fig. q, generating correlations between SNi(k,t) 
and 6Ni(—k,t). Mathematically, the consequence for 
the stabilityof the homogeneous solution is that the two 
equations ( JI9| ) effectively decouple and that only the FH 
becomes unstable at the threshold. 
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B. Correlations below shot- noise 



If the normal ordered variance becomes negative 



Once correlations between intensity fluctuations are 
identified, it is interesting to investigate if they are con- 
nected to nonclassical states of the intra-cavity fields. A 
coherent field obeys Poissonian photon statistics, which 
implies that the variance and the mean of the pho- 
ton number operator N are equal. Let us consider the 
photon number operators associated with the sum and 
difference of the intensities at different far-field points 
Ni(k) ± Nj{k'), where Ni(k) = a\(k)ai{k) and Nj(k) = 
cij(k)dj(k) are the number operators of two states a,i(k, t) 
and a,j(k, t). Since we will consider equal time correlation 
functions in the steady state of the system, from now on 
we will drop the time argument of the field operators. 
Taking out the special case i = j and k' — k which will 
be treated separately, the variance expressed in normal 
order (indicated by dots) reads 

Var^fc) ± Nj(k')) =:Vax{Ni{k) ± Nj(k')) : 
+ (:N i (k):)[a i (k),al(k)} 
+ (;N j (k'):)[a j (k'),a](k% (27) 

where Var(X) = (X 2 ) — (X) 2 . For a coherent state, the 
normal ordered variance vanishes, and the mean, given 
by last two terms in (p7j), represents the shot noise level 
for the considered quantity 



:V&i(N i ±N j ):< 0, 



(29) 



S.N. = {:Ni (A) :) [a. t (k), a] (k)} +{:Nj (k') :) [a, (k'), a] (k')] . 

(28) 



the variance becomes less than the mean, indicating sub- 
Poissonian behavior. Such a nonclassical state is iden- 
tified when the correlation normalized to the shot-noise 
level, defined as 



C£>(k,k') 



:Vax[Ni{k) ± Nj(k')]: 



: Ni(k) :) [at (k), &\ (k)] +(:N j (k') :) [a 3 {k% a) (k')} 



1, (30) 



is such that C\a (k,k') < 1. The computation of this 
quantity requires to write the normal ordered quantities 
appearing in (|30| ) in terms of anti-normal ordered quan- 
tities, since these are the quantities that are computed as 
averages in our Langevin equations associated with the 
Q representation. Using the identities 

:Ni{k): = a t (k)a\{k) =:Ni(k): +[a l (k),a\(k)}, (31a) 

:Nf(k)\ = ai(k)ai(k)at(k)al(k) 

= :N l {kf:+A:N t {k): [a t (k), a\(k)} 

+2[a t {k),a\{k)] 2 (31b) 

with three dots indicating anti-normal ordering, Eq. ( |30| ) 
reads, when expressed in terms of anti-normal ordered 
quantities 



C\f\k,k') 



'■Nax[N i {k)±N j {k')] 



CNi(k):)[ai(k),aKk)} - (:N j (k r fi[a j (V),ti(k r )] 



CNiikftiaiWraUk)] + (:N 3 (k')l)[a 3 (k'),aUk')} - [a^k) , a\{k)] 2 - a]{k')f 



(32) 



Then e.g. the normalized correlation Var[iVi(fc) ± 
Ni(— k)]/S.N. may be found by setting i = j = 1 and 
k' = —k. Equation j32|) is valid for k,k' ^ 0, while the 
special case k — will be addressed in the specific cases. 



V. LINEARIZED CALCULATIONS BELOW 
THRESHOLD 

Below threshold, the linear approximation scheme al- 
lows to derive semi-analytical expressions for the corre- 
lation functions defined in the previous section. These 
may be expressed in terms of the auxiliary correlation 
function 



where the superscript Q indicates that the average is done 
with the Q representation, corresponding to anti-normal 
ordered quantities, as indicated in the first line of ([53]). 

The starting point of ou r analysis is the set of linearized 
Langevin equations (23a) which have the exact solutions 



/ Pi(k,t) 

p 2 (k,t) 



= e M(fe)t 



/ /MM) 
/?!(-£, o) 

( m(k,t') 



M(fc)t / dt i e -M(k)t' 



\^r)* 2 {-k,t') 



(34) 



Cg(M') = (^Niik^SNjik',^), i,j = 1,2 
= {\Mk,t)\ 2 \A 3 (k',t)\ 2 ) 

-(\Mk,t)\ 2 )(\Aj(k',t)\ 2 ), (33) 



The first term in Eq. ( p4| ) describes how the intra-cavity 
fields with arbitrary initial conditions relax to the steady 
state solution and it does not contribute to the steady 
state correlations. The second term in Eq. (E3) gives the 
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response of the intra-cavity fields to the vacuum fluctua- 
tions entering the cavity through the partially transpar- 
ent input mirror. Starting from Eq. (|34|), it is possible 
to derive semi-analytical expressions for the correlations 

C?.(k,k') = {\(3 i (k,t)\ 2 \p j (k',t)\ 2 ) 

-mk,t)\ 2 m(k',t)\ 2 ) 

+ 2Re{A*A j ((3 i (k,t)/3*(k',t)) 

+AtA*(p i (k,t)/3 j (k',t))}S(k)5(k'), (35) 

where Re{-} denotes the real part. Whereas the first two 
terms in the right-hand side (r.h.s.) of Eq. ( |35| ) measure 
the correlations in the intensities of the fluctuations, the 
last two terms can be traced back to interferences be- 
tween the fluctuations and the homogeneous component 
of each field. Since these interferences only contribute 
to the equal time correlations when k — k' = 0, we will 
first concentrate on k, k' ^ and come back later to this 
special case. Henceforth, unless otherwise specified we 
consider the case k, k' ^ 0. 

The Gaussian character of the fluctuations in this lin- 
earized Langevin model allows to factorize ( |35| ) in terms 
of second order moments of the field fluctuations 

Cg(fc,fcO = |(A(fc,t)/3;(fc',t)>| 2 + |(/3 4 (fc,i)/3,(fc',t)>| 2 . (36) 

The field correlations (I3i(k,t)j3* (k' , t)) and 
(Pi(k,t)Pj(k\t)} can be best evaluated for the solu- 
tion Eq. ( |34| ) if we introduce the set of eigenvectors 
{v^) (fc)};=i 4 of the matrix M(fc), defined through: 



An arbitrary 4-component vector w can be decomposed 
on this basis 



wi(k) 
Wi(k) 



£V')(fc)vW(fc), 



(38) 



i=i 



and its components in the new basis are calculated 
via the linear transformation 



,(0 



(k) = 22 Ti m (k)w m (k). 



(39) 



This involves a 4 x 4 matrix T; m (fc) calculated as T(fc) = 
V(fc) -1 with Vi m (k) — v\ m \k). Decomposing now the 
noise vector appearing on the r.h.s. of Eq. ( |34| ) on this 
basis 



m(k,t') 

Vt(-k,t') 
allows to rewrite Eq. 



/ Pi(k,t) 
0t(-k,t) 

V/3 2 *(-M) 



= Y / V {l) (k,t')v^(k), (40) 
1=1 

41) in the large time limit as 



— / dt' 

nth Jq 



x ^ e A< ! >M(^(0 (M ' )v (0 (fc) . (41) 



1=1 



(37) The needed field correlations are given as 



M(fc)v^(fc) = A (() (fc)vW(fc). 

(f3 i (k,t)(3*(k',t)) = — ( dt' ( dt" J2 vg-iiQvftl^kS^W*^ (42a) 

n th Jq Jo ,i 

l,m= 1 

dt' ( dt" 4^ 1 (fc)4 , ;! ) 1 (fc')e A(!)(fc)( ^ t ' ) e A(m>(fc ' )(t "*'' ) (?7 (i) (fc^')^ (m) (fc'^'')>- (42b) 



(f3 t (k,t)^(k',t)) 



n t h Jo 



Lm=l 



The noise correlations in the new basis r — t t — —Lt t m t t 

fa<'>(fc,f)»7 (m) *(*',0> and (nW(k,t')n^(k',t")) are - Ui±^ 2 HUmi + mimi 

{v «Hk,t')v^*(k',t")) = A lm (k)S(k-k')S(t'-n (43a) ~^T m2 + 7 T i3 T m4 + 7 T M T m3 . (44b) 

(n {l) {k,t')n {m) {k',t")) = B lm (k)6(k + k')S(t' - t"). (43b) Inserting Eqs. @ in Eqs. (||) we can easily carry out 

the time integration, and neglecting transient contribu- 

where the matrix elements of the 4x4 matrices A(fc) tions, we end up with the following expressions 
and B(fc) can easily be evaluated in terms of the matrix 

elements T lm = T lm (k) as lim t )fi( k ' t )\ = J-Gt\k)5(k - k'), (45a) 

t^oo 3 n t h 

A lm {k) = T n T* ml - 4^T U T* 2 + Ti2V n2 lim t)) = —G ( +\k)8(k + k'). (45b) 

l t=roo n t h 

A* 

— Y T n T mi + ! T i3 T *n3 + T^T^, (44a) with 



9 



Z=l m=l 
4 4 



-(AW(fc) + AM 



4 } -i(fc)^7-i( fc ) 

-(AW(fc) + AW(fc))' 



Z=l m=l 



(46b) 



In terms of '(fc) and Gj+ } (fc), Eq. @ is given by 



C?Ak 1 k') = ^ r [\G ( i7\k)\H*(k-k') 



Ui 



+\G%\k)\ 2 5 2 (k + k>)]. 



(47) 



A. Intensity fluctuation correlations 

It is now easy to compute the normalized correlation 
function Eq. (|2q). This involves taking into account the 
commutation relation Eq. (fil) which reads 



1 



[Ai(k,t),A](k',t)] = 5 l3 —6(k - k'), 



(48) 



after rescaling space and time according to Eq. (|l2|) and 
the operators similar to the c-number fields in Eq. (14). 
We finally find 



\G ( ^(k)\ 2 8{k + k'f 



+ 



(49) 



with 7 h {k) = G ( 3j ] {k){G { 3j ] {k) - 1/2), the -1/2 in the 
parenthesis reflecting the conversion from anti-normal 
to direct ordering. Unlike the mathematical expression 
( [49| ) derived for an ideally infinite system, the correla- 
tion functions determined from the simulations will have 
peaks of a finite width, which will be determined by the 
discretization in fc-space used in the numerical codes, i.e. 
the inverse of the total length of the system. This differ- 
ence, however, will not alter the only relevant informa- 
tion, which is the height of each of these peaks. In fact, 
the quantities 



Gjj(k, —k) 



(50a) 
(50b) 



characterize the strength of the correlations between the 
modes [w](k) and [oj](-k), [2o>](ft) and [2co](-k), [u](k) 



Cg\k,uk) = 



and \2ui}(k), and [to](k) and [2w](— k) respectively. One 
easily checks that C^(k, k) = 1, as a result of an auto- 
correlation. 

All the expressions derived so far are only valid for 
nonvanishing transverse wave numbers. At k = k' = 0, 
we already observed that there are extra contributions 
to the equal time correlation function, as expressed by 
Eq. (|35|). Furthermore, in the framework of an expan- 
sion in the small parameter it is obvious that 
these extra terms even dominate, since they scale with 
\f3i(k,t)\ 2 ~ 2 1 nth, whereas the contributions on the 
first line of Eq. Q scales with \f3 l {k,t)\ 4 ~ (2jn th ) 2 . 



Hence, in the leading order, the correlation function at 
k = k' = is given by 



^12 | fe=fc/=o 



2,5(0) 
nth 



2Re(^^ 2 G^ 2 ) (0) 



-AlAZG£\oj)6(k)6(k')\. 



lfc=fc'=0" 



(51) 



Similar calculations as before allow us to derive the fol- 
lowing expression for the value of the normalized cross- 
correlation at k = k' = 



Re ( A* 1 A 2 G{~ ) (0) + AtAZGg* (0) 



(+)/ 



. (52) 



where = IAf(G^(0) - 1/4) + Re{A*< G)V(0)}- 



B. Nonclassical photon number variances 



The photon number variances considered in Sec. IV B 



can be calculated in terms of the auxiliary functions 



Gtf\k) and G\J'(k) as well. The anti-normal ordered 



•»(+) 



quantities in Eq. ( p2[ ) can be directly calculated by av- 
erages in the Langevin equation, so below threshold the 
anti-normal ordered variance is for k =/= 

;Var[A, ; (fc) ± N 3 (-k)}\ = Var[|&(fc, t)| 2 ± \^(-k,t)\ 2 }. 

(53) 

Using the commutation relations ((48|), the commutators 



in Eq. (32) are [dj(k), dj(k)] = 5(0)/n t h, and the normal- 
ized self-correlations takes the form 



C ( j f ) (k ) -k) = 



2(|G(T) (fc )| 2± |GW(fc)| 2 )-G(T)( fc ) 



G { rMk)- 1/2 



Similarly, the cross-correlations are 



G 



^(k)\ 2 ±2\G[r\k)\ 2 



T.^ik) 



^G^(k) 



v = +1,-1 



(54) 



(55) 



10 



When k 



k' 



Eq. (B3|) is no longer valid. Instead, following the procedure outlined for the normalized 



correlations we have to the leading order 0(n 



th 



Cg°(0,0) = 4- 



Re 



V A* 2 G (+) 

^3 ^3 33 



{<S)± 2A\ (a* 2 G { +\o)-\ 



(-) 



A2G 12 



(0) 



1. 



(56) 



The self-correlations become 

c 1 ( 7 ) (o,o) = o, 



-i2^ Ai/ nr(+) 



*<BT'(0) 



(57a) 

+ 4G^ } (0)-1. (57b) 

is the phase of Aj. Note that (7^(0,0) 

is actually Var [jV, (0)1 normalized to shot-noise. The 
result of Eq. (57a) is simply because the correlation 



Cj+ 3 (0,0) =4Re 



where 



Cjj \k,k') amounts to calculating the variance of zero 
for k = k' = 0. 



u 





-3 -2 -1 



1 2 3 



-3 -2 -1 



1 2 3 



u 




-0.2 
-0.4 




-0.2 
-0.4 



-1 



2 3 -3-2-10123 

k k 
FIG. 7. The linear self-correlations a) Cn(&, — k), b) 
C22(k,—k) and linear cross-correlations c) C" 2 (k,k), d) 
C" 2 (k,—k) as function of the transverse wave number for 
E/E t — 0.99. The points are numerical results while the 
lines are analytical results. 



VI. CORRELATIONS BELOW THRESHOLD 

The linearized results of Sec. gives an analytical in- 
sight to the behaviour below threshold for pattern for- 
mation. However, very close to the threshold this linear 
approximation breaks down because of critical nonlinear 
fluctuations, and additional contributions may emerge as 
for example shown in a vector Kerr model by Hoyuelos 
et al. H^] . Such nonlinear correlations can be calculated 
through numerical simulations of the full nonlinear evo- 
lution equations. 

In this section we present numerical results obtained 
from simulations of the nonlinear equations ( |l5| ) below 



threshold, with the parameters discussed in Sec. HI. Our 
numerical results are compared with the analytical re- 
sults of the previous section, and therefore also serves as 
a cross-check of our analytical and numerical methods. 



A. Linear correlations: Analytical and numerical 
results 



We first consider the strength of the correlations be- 
tween symmetric points in the far fields below the thresh- 
old for pattern formation. In Fig. [?], the four quantities 
defined by Eq. (^||) are plotted. The data are obtained 
from numerical simulations and from the analytical re- 
sults of Eqs. (|50| ) and (|52|). Very good agreement is found 
between numerics and analytical results. 

There are three main features to be considered in the 
results of Fig. |7|. First, all curves present a distinctly 
peaked behavior around the critical wave number k c for 
pattern formation, which means that the corresponding 
modes are more strongly correlated than the modes at 
any other wave number. Manifestly, this behavior is 
connected with the pattern formation mechanism and is 
closely related to the phenomenon of quantum images 
Jl0| ]. Secondly, we also note that in all four plots the 
correlations show a jump at k = 0. In a) and b) it is 
the trivial manifestation of an auto-correlation, since for 
k = 0, k and — k coincide, while in c) and d) the jump 
is due to the extra interferences with the homogeneous 
background fields as predicted from Eq. ([32]) . Finally, we 
observe that the peaks localized around k c are superim- 
posed onto smooth correlation profiles. 

The strong correlations appearing between the modes 
associated with wave numbers around k c indicate a 
strongly synchronized emission of photons in the modes 
[u](+k), H(-fc) and [2u)](+k), [2u](-k). This be- 
havior reflects the direction of instability of the sys- 
tem. As a matter-of-fact, regardless that all trans- 
verse modes of both fields are equally excited by the 
vacuum fluctuations entering the cavity, the fluctua- 
tions of the intra-cavity field modes around the criti- 
cal wave vector will be less damped than the fluctu- 
ations in the other modes. The closer to the thresh- 
old, the more the behavior of the intra-cavity fields 
will be dominated by the mode that becomes unsta- 
ble at the threshold and gives rise to the pattern. In 
the 4-dimensional phase space spanned by the fluctu- 
ation amplitudes {/3i(fc, t), f3l(-k, t),(3 2 {k, t), /3| (—k, t)}, 
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this mode is characterized by a vector with a given direc- 
tion. What we learn from the correlation functions is that 
the emerging instability results in an almost perfectly 
synchronized emission of photons in the modes [w](+fc), 
[u](-k) and [2u](+k), [2u](-k). 




E/E t 

FIG. 8. The self-correlations Cii(k c ,— k c ) (full, squares) 
and C22{kc,—kc) (dashed, circles) and the cross-correlations 
Ci2(k c , —k c ) (dotted, triangles) and C]^(fe c , k c ) (dash-dotted, 
diamonds) as function of the pump normalized to the thresh- 
old. The points are numerical results while the lines are ana- 
lytical results. 



The dominance of this particular mode when the 
threshold is approached is confirmed by the study of the 
strength of these correlations as a function of the pump. 
In Fig. ^ we follow the height of the peaks at k = k c of the 
four linear correlations displayed in Fig. f^, as a function 
of the pump level E/Et. The most immediate obser- 
vation is that all the correlations become perfect in the 
limit E — > E t . This asymptotic behavior can be under- 
stood from the linearized fluctuation analysis presented 
in Sec. fv| It is enough to observe that Eqs. (fl6|) involve 
the inverse of th e real part of the eigenvalues of the lin- 
ear system ( |23a|) . The dominance at the th resho ld of the 
undamped eigenmode of the linear system ( |23a| ) emerges 
from the fact that here the real part of the associated 
eigenvalue precisely goes to zero. Thus, the decrease in 
the correlations as we move away from threshold can be 
seen as the result of the coexistence of different eigen- 
modes. Physically the emergence of these correlations is 
much less intuitive than the ones in an OPO. As a matter- 
of-fact, in the OPO below the threshold momentum con- 
servation is enough to predict the existence of correla- 
tions between the fluctuations in the modes [o;](+fe) and 
[ui](— k). In the presence of the 4-mode interaction of 
SHG, the momentum conservation gives a global condi- 
tion involving all four beams (at [w] (+£;), [w](— k) and 
[2w](+fe), [2u](—k)). These correlations in fact arise in 
connection with the emergence of an instability. 

Turning now to the cross-correlation between the ho- 
mogeneous components of the fields, we observe that 



Ci 2 {k — 0, k' = 0) in Fig. [7] is negative, reflecting an 
anticorrelation of the photons associated with the FH 
and SH homogeneous waves. In other words, the cre- 
ation of a photon [2u>] (0) implies the destruction of (two) 
photons [w](0) and vice versa. The origin of this corre- 
lation is much simpler to understand than the previous 
one: The two modes [w](0) and [2u>] (0) being macroscop- 
ically populated, the vacuum fluctuations simply induce 
transitions between these two modes, according to step 
1) in the scheme in Fig. ||. In Fig. || we plot this correla- 
tion as a function of the pump. Comparing the value of 
the correlations below and above threshold, we observe 
that very close to, but below, the threshold, the tendency 
of the curve is reversed and it anticipates the behavior 
of the correlation above threshold. These are nonlinear 



correlation effects that will be discussed in Sec. VI B 




FIG. 9. The linear cross-correlation Ci 2 (k = 0, k' = 0) as 
function of the pump normalized to the threshold, comparing 
numerical results (points) with the analytical result (line). 
Open (full) symbols are numerics below (above) Et- 



Finally, we would like to discuss the smooth contribu- 
tions to the correlations displayed in Fig. fj]. We first note 
that these are not connected with the pattern instability. 
This was checked by considering very low pump values 
for which the peaks around k c completely vanish, while 
the smooth structures of the curves remain. Consider- 
ing the central region of the curves, roughly for \k\ < k c , 
the most striking observation is the absence of correla- 
tions between the fluctuations in the modes [w](fe) and 
[2w](fc), whereas [oj](k) and [2w](— k) are correlated, as 
well as [2w](A:) with [2w](— k). This behavior seems to 
indicate the existence of a symmetry restoring principle 
in the dynamics of the intra-cavity fields. As a matter- 
of-fact, the absence of correlations between [w](/c) and 
[2w](fc) implies that the fluctuations of the numbers of 
pair productions through step 2) and the fluctuations 
of the number of conversions — > [2w](fe) through 
step 3) occur independently of each other. However, 
while step 2) of Fig. || conserves the k — > —k symme- 
try of the system, step 3) does not. As a consequence, 
a positive fluctuation in the number of times step 3) oc- 
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curs ([u;](fe) + [w](0) — > [2w] (fc)), automatically implies 
that there will be less [w](fc) than [oj}(— k) in the system, 
and more [2w](fc) than [2w](— k). The correlations ob- 
served may indicate that the system will try to restore 
the k — > — k symmetry by down-converting [2cj](0) — > 
[cj](fc) + k), producing a surplus of [w](— k) which 
again will produce more [2w](— k). These mechanisms 
seems to fit well with the relative strengths of the cor- 
relations observed in the central region of Fig. g. The 
strongest is always Cj\(fc, — k), in agreement with the 
fact that the twin photon emission is the principal source 
of correlations in the system. Weaker is the correlation 
Cf 2 (k, —k) and even weaker C% 2 (k, —k). This interpreta- 
tion is consistent with the way the correlations at k = k c 
depart from the value 1 at threshold, when the pump is 
lowered, as displayed in Fig. ||[ 



level is taken beyond threshold, cf. Sec. VII. Therefore 
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FIG. 10. Photon number variances for E/Et = 0.99 show- 
ing (k, — k) (full, diamonds) and (k, — k) (dashed, 
squares). The lines are analytical results while the points are 
numerical simulations. The shot-noise level C = 1 is indicated 
by a thin dotted line. 



We now turn our attention to the study of the fluctu- 
ations in the sum and difference of the photon numbers 
at symmetrical points of the far field. We first consider 
the twin beam photon variances for the FH, C^'{k, —k) 
defined in Eq. (|30| ) and shown in Fig. [HJ The results 
are symmetric with respect to the substitution k — > —k, 
wherefore we plotted this quantity for positive k, shifting 
the origin for better view of the specific behavior at k = 0. 
The linearized calculation predicts sub shot-noise statis- 
tics in the difference N\(k) — Ni(—k) for all wave num- 
bers. For large wave numbers the analytical result for the 
correlation approaches the value 1/2. It is interesting to 
keep in mind that for the OPO, the same quantity is equal 
to 1/2 independently of the wave number |3^,Q. In the 
SHG case additional processes taking place in the cavity 
result in a smooth fc-dependence of c[i\k, —k). These 
characteristics do not depend much on the value of pump, 
and are not changed significantly even when the pump 



the statistics of the intensity difference is not directly 
affected by the pattern formation mechanism. A rad- 
ically different situation occurs for the sum-correlation 
(k, —k), which shows a strong peak around k = k c . 
For the pump value used in Fig. |l0] the peaks correspond 
to a maximum value C[i\k c , —k c ) ~ 35. This behav- 
ior is connected with the increase of the fluctuations in 
the modes associated with the pattern instability when 
the threshold is approached, leading to a large excess 
noise in the statistics of the intensity of the individual 
modes [uj]{k) and k). This excess noise in each in- 
tensity cancels when the difference N\(k) — N\(—k) is 
considered leading to a sub-Poissonian statistics, while 
it is still present in the sum Nx(k) + Nx(— k). For large 
k the correlation approaches 1.5, coinciding again with 
the corresponding value for the OPO. Finally, as before, 
the jumps at k — are due to contributions from the 
homogeneous steady states, cf. Eqs. (|56|)-(|57|). 
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FIG. 11. Photon number variances for E/Et = 0.99 show- 
ing {k, — k) (full, diamonds) and C^\k,—k) (dashed, 
squares) . 



The corresponding photon number variances 
C^\k,-k) for the SH field are shown in Fig. [ll]. In con- 
trast to the FH correlations there is almost no sub-shot 
noise behavior in the difference correlation (k, — k). 
In other words, the SH beams only display very weak 
nonclassical correlations. As for the FH field, the 
emerging instability does not influence the noise level 
in C^\k, —k), but C^{k, —k) displays a large amount 
of excess noise in the vicinity of k c . The asymptotic 
large k behavior for both correlations (k, — k) and 
C^{k, — k) is analytically found to correspond to the 
shot-noise limit 1.0. 

The cross-correlations c[ 2 \k,k) and cj^(A;,A:) are 
shown in Fig. |l2|. The linearization approach predicts 
that these correlations are always above the shot-noise 
limit. Furthermore, at small wave numbers we note that 
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the variances of the sum and difference coincide. This can 
only occur when the fluctuations in the individual modes 
[u>](k) and [2wl(fc) are uncorrelated, what was indeed ob- 
served in Fig. 0. Moreover, both the sum and difference 
correlations show a large excess noise at k = k c , which is 
slightly weaker for the difference, as the result of a partial 
noise cancellation. 
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cal situation is related to what was observed in the nor- 
malized correlations where the cross-correlation between 
N%(k) and N 2 (~k) was stronger than the almost van- 
ishing cross-correlation between Ni{k) and N 2 (k). At 
k = k c a large amount of excess noise dominates the be- 
havior of both the sum and the difference correlation and 
the two correlations show a pronounced peak. For large 
k the correlations approach the shot noise limit, as seen 
for the other cross-correlations in Fig. [12[ 

Olsen et al. |3^] have investigated the system without 
spatial coupling corresponding to our results at k = 0, 
and they find that, for certain detunings, the variance of 
the sum of the FH and SH intensities are more strongly 
quantum correlated than the variance of the individ- 
ual intensities, due to the anti-correlation between them. 
Var[iVi(0)]/5.AT. and Vax[N 2 (0)]/S.N. can be seen from 
Fig. [h] and O, respectively, at k = 0. Both are larger 
than the Var[iVi(0) + N 2 {0)]/S.N. observed in Figs. |§ 
and H3, so that our results confirm the ones of Ip5[] . 



B. Nonlinear correlations: Numerical results 



FIG. 12. Photon number variances for E/Et = 0.99 
showing C\^{k, k) (full, diamonds) and C^\k,k) (dashed, 
squares) . 
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FIG. 13. Photon number variances for E/Et = 0.99 show- 
ing c'^(fc, — k) (full, diamonds) and Cj^ (k,—k) (dashed, 
squares) . 



The cross-correlations C\ 2 \k,—k) and 



C[ 2 \k, - 



■k) 



are shown in Fig. |12|, and here the difference correlations 
interestingly go below the shot-noise limit as long as k 
is not too close to the critical wave number. It is worth 
pointing out that the difference Nx(k) — N 2 (—k) shows 
nonclassical behavior while the difference Ni(k) — N 2 {k) 
(shown in Fig. Il2]) does not. This somehow paradoxi- 



So far we have only considered the correlations pre- 
dicted by the linearized equations. In order to go be- 
yond this regime, we use our numerical simulations to 
search for nonlinear fingerprints in the correlations and 
in particular for the emergence of new correlations, i.e. 
Cfj(k, k') with k ±k'. Of particular interest is to look 
for correlations between the homogeneous steady states 
(k = 0) and the states with k — ±fc c , Gy (0, ±k c ). From a 
technical point of view this task turned out to be difficult 
because nonlinear contributions to the correlation func- 
tions were only observable for pump values extremely 
close to threshold, in a region where the characteristic 
time of the dynamics diverges because of critical slowing 
down. This translates into very long transients and the 
need of equally long simulations. 

We have observed some indication of nonlinear corre- 
lations for a pump E/Et = 0.99999, which became very 
clear when using E/E t = 0.999999. For this value of the 
pump, we show in Fig. 14 our results for C" 2 (fc, kl = 0) 
and C" 2 (fc' = 0, k): these curves put into evidence an 
anti-correlation between the modes [w](±/c c ) and [2u](0), 
and between [2cu] (±fc c ) and [u>] (0) . They present a very 
sharp peak structure, with a width determined by the 
distance between two adjacent points of the discretized 
fc-space used for the simulations. This is due to the fact 
that we now consider the correlation functions at fixed k 1 
and let k vary. These correlations arc a result of nonlinear 
amplification of the diverging fluctuations as the thresh- 
old is approached. The negative nature of the correlation 
is connected with the fact that the fields with nonzero 
average values (here the homogeneous components) act 
as a " reservoir" of photons for all processes occurring in 
the cavity. As we will show later, they are a precursor 
of the behavior of the correlations above the threshold. 



14 



The correlations at k = correspond to the linear cor- 
relation shown in Fig. ||. The bottom plot in Fig. |l4| 
shows the nonlinear correlations C 4 "(fc = +k c ,k' = 0) 
as the threshold is approached. The correlations are 
nonzero only for E / E t > 0.9999. As the nonlinear cor- 
relations set in the nonlinear channels in step 2) and 3) 
of Fig. |^ become stronger and this weakens the correla- 
tions induced by the channel of step 1), which is exactly 
what we observed in Fig. ^; Cf 2 (0, 0) becomes less corre- 
lated very close to the threshold. Moreover, we see that 
the correlations Cj\(0, +k c ) and C" 2 (0, +k c ) have almost 
identical values, and the same holds for C 22 (+fc c ,0) and 
C™ 2 (+k c , 0) . This interesting behavior can be traced back 
to the fact that close to the threshold the fluctuations 
8Ii(k c ) and <5/ 2 (fc c ) are perfectly correlated, as displayed 
by Fig. whereas the slight anti-correlation between 
<5/i(0) and <5/ 2 (0) is responsible for the lower values of 
C 2 l 2 (+fc c ,0) and Cf 2 (+fc c ,0) with respect to C{\(0,+k c ) 
andC^O.+fcc). 
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FIG. 14. Above: Nonlinear cross-correlations 

C" 2 (k,k' = 0) (left) and Of 2 (7c' = 0,k) (right) as function 
of k for E/E t = 0.999999. Below: Semi-log plot the nonlin- 
ear correlations Cj™ (k — +k c , k' — 0) as function of E/Et. 



VII. CORRELATIONS ABOVE THRESHOLD 

Above the threshold for pattern formation the lin- 
earized equations ([l9|) are no longer valid. As displayed 
in Fig. ^, above the threshold not only the homoge- 
neous modes, but also all modes with wave numbers 



k = ±fc c , ±2fc c , ±3fc c , . . ., will present a macroscopic pho- 
ton number. Linearizing around the steady state pat- 
tern solution above the threshold under the assump- 
tion of small fluctuations, one obtains new linear equa- 
tions for the far field fluctuation amplitudes, which take 
into account three-wave processes such as [2w](fc c ) <-» 
[u)](k) + [u)](k c -k) or [2u](k) <-> [u](k c ) + [w](k-k c ). In 
analogy to the situation below the threshold a linear fluc- 
tuation analysis above the threshold predicts, in addition 
to the correlations already present below the threshold, 
the existence of additional correlations between the fluc- 
tuations 5I\(k) and 8I\(k c — k), and between <5/ 2 (A:) and 
81\{k — k c ). We will not report here the explicit results 
of this cumbersome linear analysis and restore directly 
to the numerical analysis of the full nonlinear Langevin 
equations. 
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FIG. 15. The self-correlations a) C^(k, ~k), b) CJ 2 (fc, -k) 
and cross-correlations c) Ci 2 (k,k), d) Cx 2 (k, — k) as function 
of the transverse wave number for E/Et = 1.05. 



To investigate the implications of the new field con- 
figuration above the threshold on the intensity correla- 
tions, we first consider the correlations C™ (fc,fc'). The 
same normalized correlations discussed in Fig. ffl below 
the threshold are plotted in Fig. [l5| for a pump value 
above the threshold. We observe that the correlations at 
k = ±fc c decrease from their threshold value and are no 
longer perfect as they were at the threshold. A closer 
look actually reveals a dip in the correlations exactly at 
the pixels corresponding to k = ±k c . A tentative expla- 
nation for this is based on the fact that now the modes at 
the critical wave number have a finite average value, con- 
nected with macroscopic photon numbers in these modes, 
whereas the neighbouring pixels are significantly less pop- 
ulated, cf. the far field of Fig. || In comparison the 
normalized correlations C 1 ™- (fc, k') show a much smoother 
behaviour around k r . Hence, the observed reductions in 
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the correlations above threshold at k = ±k c are con- 
nected with spontaneous population exchanges between 
these macroscopically populated modes. 

In Fig. H the peaks at k = ±k c of Fig. [l5| are fol- 
lowed as function of the pump. The behavior is very 
similar to what is seen below the threshold. Close to the 
threshold the correlations are perfect, and as the pump 
is taken further away from E t the correlations become 
weaker. Below the threshold this was explained through 
an eigenvalue competition, while above the threshold the 
explanation is that the competitions between the states 
become stronger. 

The k = cross-correlation is plotted in Fig. |9|, and 
above the threshold there is a loss of anticorrelation or 
there is even a small positive correlation. This might be 
attributed to the macroscopic and independent occur- 
rences of the processes of step 2) and 3) in Fig. |[ 
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FIG. 16. The correlations C"j(k = +k c , k! — 0) as function 
of the pump relative to the threshold. The gray symbols are 
the correlations below the threshold from Fig. h4. 







1 1 1 1 1 1 1 1 1 1 


1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 

□ 

□ 




□ 

□ 




□ 

n □ 

c tt, □ °n : 

^td n □□ 


□ 




^ooooooo ooooooo** 5 


oo°° °o„ 

oooooooo* 


' .... 1 .... 1 





-0.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 
k 

FIG. 17. Photon number variances for E/Et = 1.05 show- 
ing C[i\k, — k) (diamonds) and C[~^\k, —k) (squares) from 
a numerical simulation. 



We saw in Sec. VI B nonlinear correlations just be- 
low the threshold, and in Fig. [If] the peaks correspond- 
ing to these correlations are plotted in order to follow 
the progress above the threshold. The strongest anti- 
correlation is observed just above the threshold, E/E t — 
1.0001 and as the pump is increased the correlations be- 
come weaker due to increasing competition of processes 
involving higher harmonics. Moreover, the connection 
between the self-correlations and cross-correlations seen 
below E t only remains very close to the threshold, so 
as the pump is increased C{\(0, +k c ) ^ C™ 2 (0, +k c ) and 
C2 2 (+k c ,0) ^ C[ l 2 (+k c ,0). This is related to the loss of 
perfect correlations away from the threshold. 

In Fig. [I?] the photon number variances c[f\k,— k) 
above the threshold are presented. Comparing these re- 
sults with the corresponding ones below the threshold 
from Fig. [n] we observe that they are very similar. Gen- 
erally, the correlation c[i \k, —k) does not change much 
with the pump level, and this fact has also been observed 
in the OPO @. The sum correlation c[p(k, -k), how- 
ever, contains peaks that are very sensitive to the pump 
level, both below and above the threshold. The behavior 
discussed here for the FH is also valid for the SH and the 
cross-correlations. 



VIII. CONCLUSION AND DISCUSSION 



We have used the master equation approach to de- 
scribe the spatiotemporal dynamics of the boson intra- 
cavity operators in second-harmonic generation, and we 
included in the model quantum noise as well as diffrac- 
tion. Our study is based on the Q representation to de- 
scribe the dynamics of the quantum fields in terms of a 
set of nonlinear stochastic Langevin equations for equiv- 
alent c-number fields. The choice of the Q representa- 
tions gives some restraints on the parameter space, but 
we have checked that similar results are obtained using 
the approximated Wigner representation in other param- 
eter regions. 

A simple scheme describing the microscopic photon 
interaction that underlies the process of pattern forma- 
tion has guided us in our analytical and numerical stud- 
ies of the spatial correlations. Equal time correlations 
between intensity fluctuations were used to investigate 
the strength of the correlations between different modes. 
Also, possible nonclassical effects, such as twin beam 
correlations, were considered by calculating the photon 
number variances of the intensity sums and differences 
between spatial modes of the FH and SH fields. 

We have found that at the threshold for pattern for- 
mation the Fourier modes with the critical wave number 
are perfectly correlated for the FH field, the SH field and 
also between the FH and the SH field. As the distance 
to the threshold is increased these correlations become 
weaker, which was shown analytically to be due to the 
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competition of the eigenvalues of the linear system de- 
scribing the system below the threshold. At large wave 
numbers, only the correlation between opposite points 
of the FH far field survives. This correlation is always 
found to be stronger than the others, which is consistent 
with the fact that the twin photon emission at the funda- 
mental frequency is the primary source for correlations 
in the system. For far field modes around the critical 
wave number the self-correlations as well as the cross- 
correlations between FH and SH photons are linked to 
the pattern forming instability. 

Very close to the threshold the linear analysis breaks 
down. The numerical simulations below the threshold 
showed the existence of nonlinear correlations which in- 
volve the k — mode and these are also seen above 
the threshold. The other correlations described above 
are also found above the threshold, but their strength 
decreases when moving away from the threshold. This 
can be understood from the fact that additional pro- 
cesses come into play, mainly consisting in population 
exchanges between the macroscopic fields at the critical 
wave number and its harmonics. 

The intensity differences between opposite points of 
both the FH and SH far fields, as well as the cross- 
correlation between the two have been shown to exhibit 
nonclassical sub-shot noise behavior. These properties 
for the intensity difference turn out not to be sensitive to 
the process of pattern formation, since the correspond- 
ing correlations depend very weakly on the distance to 
the threshold and show no particular structure close to 
the critical wave number. The emerging pattern is con- 
nected with increased fluctuations in the modes with 
wave numbers around the critical wave number, leading 
to an excess noise in the corresponding individual inten- 
sities. Therefore, the sub-Poissonian statistics of the in- 
tensity differences reveal a partial noise cancellation. On 
the contrary, the sum of intensities clearly exhibit peaks 
around the critical wave number, originating from excess 
noise connected with the formation of a pattern. 

In this work we considered equal time correlations cal- 
culated for the intra-cavity fields. This approach turned 
out to be very useful to understand the intra-cavity field 
dynamics. For the output fields we expect that the non- 
classical correlations of the intra-cavity fields will remain 
below shot noise. The quantitative assessment of the 
amount of noise reduction or excess noise with respect 
to the shot noise level requires a specific additional cal- 
culation. For future work it would also be interesting to 
calculate the output fluctuation spectra at frequency 
for the difference and sum of intensities, which reflect 
the full amount of quantum correlations induced by the 
microscopic processes taking place inside the cavity, as 
for example considered for a vectorial Kerr model in [p0| . 
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